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1.0  SUMMARY 


A  primary  aspect  involved  in  the  navigation  of  a  wide  variety  of  vehicles,  such  as  ground  vehicles 
and  aircraft,  is  the  use  of  inertial  measurement  units.  The  purpose  of  the  inertial  measurement 
unit  is  to  provide  an  internal  measure  of  the  linear  acceleration  and  angular  velocity  of  the  vehicle 
in  order  to  be  used  by  a  navigation  system  in  the  dead-reckoning  of  the  vehicle’s  position  and 
attitude.  Due  to  the  nature  of  sensors,  the  measurement  of  the  linear  acceleration  and  angular 
velocity  are  corrupted  by  several  error  sources,  such  as  bias  and  noise,  amongst  others.  Another 
key  element  involved  in  the  navigation  of  vehicles  is  the  use  of  external  aiding  via  more  conventional 
sensors,  such  as  range-finders  and  star  cameras.  The  purpose  of  these  sensors  is  to  aid  navigation 
by  providing  measurements  relative  to  the  external  environment  in  which  the  vehicle  is  operating. 
In  order  to  assess  the  coupled  impact  of  these  sensors,  along  with  their  inherent  uncertainties,  on 
the  navigation  accuracy  for  an  arbitrary  vehicle,  a  simulation  and  analysis  tool  is  developed  that 
can  model  several  inertial  measurement  unit  mechanizations,  model  several  sensors  for  the  purpose 
of  external  aiding,  simulate  the  error  sources  encountered  in  the  acquisition  of  measurement  data, 
emulate  the  navigation  software,  and  perform  a  range  of  analyses  via  Monte  Carlo  simulation  and 
linear  covariance  analysis.  The  developed  tool  is  applied  to  a  ground  vehicle  navigation  scenario  in 
order  to  assess  the  performance  of  different  sensors  and  provide  conclusions  as  to  when  each  sensor 
should  be  used  to  achieve  the  best  navigation  performance. 

2.0  INTRODUCTION 

Inertial  measurement  unit  (IMU)  technology  is  used  to  drive  a  dominant  element  of  navigation 
systems  on  vehicles  ranging  from  submarines  to  ground  vehicles,  missiles,  and  spacecraft.  The 
purpose  of  the  IMU  is  to  provide  an  internal  measure  of  the  changes  in  the  linear  and  angular 
motions  of  the  vehicle  due  to  non-gravitational  stimuli  by  using  triads  of  single-axis  accelerometers 
to  measure  the  linear  motion  and  gyros  to  measure  the  angular  motion  via  measurement  of  the 
angular  velocity.  In  doing  so,  the  position,  velocity,  and  attitude  of  the  vehicle  can  be  predicted 
from  a  prescribed  starting  condition  along  with  a  representation  of  the  confidence  in  the  prediction, 
i.e.  an  associated  uncertainty. 

Traditionally,  IMU-based  navigation  is  mechanized  in  one  of  two  ways:  gimbaled/gyro-stabilized 
or  strapdown.  In  gimbaled  mechanizations,  the  triad  of  accelerometers  is  placed  on  a  “stable  table” 
at  the  center  of  set  of  interconnected  gimbals  or  at  the  center  of  a  fluidic  stabilization  system.  In 
the  strapdown  mechanization,  the  triad  of  accelerometers  and  gyros  is  fixed  (strapped  down)  to 
the  body  of  the  vehicle,  and  there  is  no  stable  table.  The  different  mechanizations  effectively  alter 
the  way  in  which  the  triad  orientation  is  handled:  for  gimbaled  systems,  mechanical  tracking  is 
used,  whereas  computational  tracking  is  used  for  strapdown  systems.  Common  to  every  system, 
the  outputs  of  the  single-axis  sensors  that  make  up  the  accelerometer  and  gyroscope  are  used  in 
conjunction  with  the  navigation  system  software  to  dead-reckon  (i.e.  using  only  the  IMU)  the  state 
(position,  velocity,  and  attitude)  of  the  vehicle. 

Another  core  element  of  vehicular  navigation  is  the  use  of  external  sensor  data,  in  conjunction 
with  the  IMU  data,  to  improve  knowledge  of  the  position,  velocity,  and  attitude  of  the  vehicle.  Data 
acquired  from  range  finders,  star  trackers,  and  line-of-sight  sensors,  amongst  others,  can  be  used 
to  provide  refinements  in  the  uncertainty  with  which  the  state  of  the  vehicle  can  be  determined. 

Using  the  IMU  to  predict  the  position,  velocity,  and  attitude  of  a  vehicle  inherently  leads  to 
uncertainty  in  the  aforementioned  localization  variables.  Even  in  the  case  where  perfect  initial 
conditions  are  known  (i.e.  there  is  no  uncertainty  in  the  position,  velocity,  and  attitude),  errors  in 
the  IMU  parameterization  will  flow  into  the  translational  and  rotational  states,  causing  a  growth 
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in  the  uncertainty  on  the  prediction  of  the  position,  velocity,  and  attitude  through  a  degradation 
of  the  associated  confidence  in  these  predictions.  Using  the  external  sensor  data  can  mitigate  this 
growth  of  uncertainty  by  adding  information  about  the  vehicle’s  state. 

Employing  external  measurement  data  increases  the  system  complexity  and  cost  by  introduc¬ 
ing  additional  system  components;  however,  if  the  system  complexity  growth  is  accompanied  by 
improved  performance  of  the  navigation  system,  then  operation  of  the  vehicle  can  be  considerably 
improved.  This  can  lead  to  more  reliable  operation  of  autonomous  or  remotely  operated  vehicles 
and  it  can  lead  to  more  robust  fault  and  failure  monitoring  for  human  operated  vehicles.  The  trade 
between  increased  complexity  and  improved  performance  needs  to  be  thoroughly  analyzed  and 
understood  in  order  to  assess  whether  or  not  a  suite  of  external  measurement  senors  is  efficacious. 

Analysis  of  IMU-driven,  externally  aided  navigation  systems  can  be  accomplished  in  several 
ways.  One  of  the  most  flexible  approaches  is  to  use  techniques  from  Monte  Carlo  simulation  and 
particle  filtering.  This  approach  allows  for  the  utilization  of  highly  complex  systems  in  a  black 
box  configuration  with  minimal  intrusion  into  the  governing  equations;  however,  this  can  be  quite 
computationally  complex  and  require  significant  computing  resources.  An  alternative  is  to  use 
linear  covariance  analysis  to  predict  the  mean  and  covariance  that  result  from  the  implementation 
of  a  certain  combination  of  an  IMU  and  other  external  sensors.  This  method  accounts  for  the  errors 
in  all  of  the  parameters,  but  does  require  linearization  to  be  performed  in  order  to  determine  the 
resulting  uncertainties.  The  benefit  of  linear  covariance  analysis  over  Monte  Carlo  analysis  is  that 
linear  covariance  analysis  is  computationally  efficient  and  can  be  performed  on  many  configurations 
of  the  internal/external  sensor  combinations  very  expediently. 

In  order  to  quantify  the  uncertainties  involved  with  a  specific  navigation  system  that  is  com¬ 
prised  of  internal  and  external  sensors,  SAIMUN  (Simulation  and  Analysis  of  IMU-based  Naviga¬ 
tion)  was  developed.  This  tool  utilizes  a  known  truth  trajectory,  IMU  parameterization  of  interest, 
prescribed  suite  of  external  aiding  sensors,  and  initial  conditions  to  calculate  performance  criteria. 
The  tool  makes  use  of  Monte  Carlo  analysis  and  linear  covariance  analysis  to  propagate  the  posi¬ 
tion,  velocity,  and  attitude  uncertainty  forward  in  time.  Monte  Carlo  analysis  is  used  to  investigate 
the  uncertainties  for  IMU-only  navigation  systems  and  to  baseline  performance  in  this  case.  Linear 
covariance  analysis  is  used  to  compare  to  Monte  Carlo  analysis  for  IMU-only  navigation  and  to 
perform  externally  aided  navigation  analysis  as  well. 

SAIMUN  was  developed  to  consider  two  different  IMU  mechanizations:  strapdown  and  space 
stabilized.  The  strapdown  IMU  mechanization  relies  on  microelectromechanical  systems  (MEMS) 
technology  for  the  accelerometer  and  gyroscope  triads.  The  triads  are  stationary  with  respect  to 
the  IMU  case;  thus,  the  attitude  is  dead-reckoned  numerically  using  the  gyroscope  data  to  provide 
the  vehicle  attitude.  The  space  stabilized  IMU  mechanization  relies  on  a  mechanical  platform  to 
track  an  inertial  reference  frame.  Additionally,  SAIMUN  was  developed  to  consider  three  types 
of  external  measurements:  range/range-rate  measurements,  stellar  line-of-sight  measurements,  and 
bearing  angles  measurements.  The  objective  of  the  tool  is  to  provide  an  analysis  of  the  feasible 
navigation  accuracy  given  some  combination  of  IMU  and  external  sensors. 

The  remainder  of  this  report  is  organized  as  follows.  Section  3.1  provides  an  overview  of  the 
simulation  environment  that  was  developed.  Section  3.2  details  the  modeling  of  the  error  sources 
which  contribute  to  the  IMU  corruption  and  provides  the  overall  models  for  the  measured  IMU 
outputs  that  are  used  in  the  strapdown  and  space  stabilized  implementations.  Following  is  a 
detailed  development  of  the  dynamical  models  used  for  the  two  IMU  configurations  in  Section  3.4. 
Section  3.5  describes  the  Monte  Carlo  analysis  technique  and  its  usefulness  in  validating  other 
uncertainty  analysis  techniques.  Section  3.6  outlines  the  three  sensors  used  for  external  aiding 
including  application  and  their  dynamical  models.  The  theory  and  method  of  implementation  of 
the  linear  covariance  analysis  is  presented  in  Section  3.7.  Section  4.0  presents  the  results  of  a  ground 
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vehicle  navigation  simulation  in  which  two  IMUs  are  considered  and  compared  using  the  developed 
tool.  Section  5.0  concludes  the  report  with  a  discussion  as  to  how  the  developed  tool  is  used  to 
select  an  appropriate  IMU  for  the  ground  vehicle  navigation  simulation  as  well  as  a  discussion  on 
other  capabilities  that  SAIMUN  offers. 

3.0  METHODS,  ASSUMPTIONS,  AND  PROCEDURES 
3.1  SAIMUN  Overview 

In  order  to  analyze  the  effect  of  inertial  measurement  unit  and  external  sensor  performance  on 
navigation  systems,  SAIMUN  was  developed.  SAIMUN  utilizes  Monte  Carlo  sampling  and  linear 
covariance  analysis  to  propagate  the  initial  position,  velocity,  and  attitude  uncertainty  forward 
in  time,  while  accounting  for  the  injection  of  uncertainty  into  the  position,  velocity,  and  attitude 
due  to  error  corruption  present  in  the  IMU.  Updates  to  the  uncertainty  are  also  considered  when 
incorporating  external  sensor  data  using  linear  covariance  analysis. 

The  effect  of  the  following  IMU  error  sources  on  the  position,  velocity,  and  attitude  uncer¬ 
tainty  is  quantified  in  SAIMUN:  startup  bias,  bias  instability,  random  walk,  scale  factor  error,  axes 
misalignment,  axes  nonorthogonality,  and  quantization  effects.  The  modeling  of  these  IMU  error 
parameters  is  outlined  in  Section  3.2.  The  location  of  the  IMU  with  respect  to  the  center  of  mass 
(CM)  and  the  IMU  timestep,  which  is  the  time  between  successive  measurements  of  the  acceleration 
and  angular  velocity  provided  by  the  IMU,  are  also  specified.  Finally,  the  mechanization  type  for 
the  IMU,  which  can  be 

•  strap  down  or 

•  space  stabilized 

is  given.  The  IMU  mechanizations  are  further  described  in  Section  3.4. 

In  addition  to  information  about  the  IMU  configuration  and  error  sources,  the  tool  also  requires 
a  tabulated  true  trajectory,  which  is  provided  as  an  input  to  SAIMUN  to  ensure  that  the  tool  can  be 
applied  to  any  vehicle,  ranging  from  ground  vehicles  to  aquatic  vehicles,  air  vehicles,  and  spacecraft. 
The  tabulated  truth  trajectory  consists  of  the  following  information: 

•  inertial  position  vector  of  the  CM  of  the  vehicle, 

•  inertial  velocity  vector  of  the  CM  of  the  vehicle, 

•  inertial  acceleration  vector  of  the  CM  of  the  vehicle, 

•  attitude  of  the  vehicle  body  frame,  expressed  as  a  quaternion  and  defined  as  a  3-2-1  rotation 
from  the  inertial  frame, 

•  angular  velocity  of  the  vehicle  body  frame  expressed  in  the  vehicle  body  frame, 

•  angular  acceleration  of  the  vehicle  body  frame  expressed  in  the  vehicle  body  frame,  and 

•  seconds  elapsed  since  the  beginning  of  the  trajectory. 

An  initial  epoch,  which  corresponds  to  the  start  of  the  tabulated  truth  trajectory  and  is  specified 
as  a  Julian  Date,  is  also  provided  to  link  the  truth  trajectory  to  a  known  reference  time. 

SAIMUN  is  also  integrated  with  the  ability  to  process  measurements  from  external  sources.  For 
each  sensor  the  following  information  is  required: 


3 

Approved  for  public  release;  distribution  is  unlimited. 


•  measurement  frequency, 

•  measurement  error  parameters  given  as  a  covariance  matrix, 

•  sensor  position  and  attitude  with  respect  to  the  IMU,  and 

•  latency  associated  with  each  measurement. 

The  analysis  and  initial  distribution  inputs  required  are: 

•  Number  of  Monte  Carlo  sample  trajectories.  This  parameter  is  a  trade  between  computational 
effort  required  and  the  accuracy  of  the  results. 

•  Covariance  of  the  initial  position,  velocity,  and  attitude  of  the  vehicle,  corresponding  to  the 
initial  epoch. 

SAIMUN  then  implements  linear  covariance  and/or  Monte  Carlo  analysis  to  propagate  the  initial 
uncertainty  forward  in  time,  while  potentially  accounting  for  the  presence  of  external  sensor  data, 
and  outputs  the  following  plots: 

•  position  standard  deviation  over  time, 

•  velocity  standard  deviation  over  time,  and 

•  attitude  standard  deviation  over  time. 


3.2  IMU  Modeling 

The  acceleration  and  angular  velocity  measured  by  an  IMU  are  corrupted  by  a  variety  of  error 
sources.  SAIMUN  accounts  for  the  most  dominant  sources  of  error  by  specifying  a  probability 
density  function  (pdf)  for  the  error  parameters  to  define  each  error  source.  The  pdfs  implemented 
in  SAIMUN  are 

•  the  Dirac  distribution,  which  is  used  to  represent  a  constant  parameter  specified  solely  by  the 
mean, 

•  the  normal  distribution,  which  represents  a  Gaussian  distributed  parameter  that  is  solely 
defined  by  a  mean  and  variance  (or  standard  deviation),  and 

•  the  uniform  distribution,  which  represents  a  parameter  that  is  equiprobable  between  lower 
and  upper  bounds. 

These  distributions  are  graphically  represented  in  Fig.  (1). 

The  IMU  model  accounting  for  error  sources  is  given  for  the  accelerometers  and  gyroscopes  as 

am,k  =  aQ,k(^{I  +  Na  +  Ma)  ( I  +  Sa )  (T™U +  bafi  +  ba ^  +  wa^)^j  (la) 

w. m,k  =  ,k({I  +  Ng  +  Mg)  ( I  +  Sg)  -g  0  _|_  ifgk  -)_  Wg:i^j^j  ,  (lb) 


where 


ang  k  true  non-gravitational  inertial  acceleration  experienced  by  the  IMU  expressed  in 

the  inertial  frame  at  tk, 
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Figure  1.  Three  Distributions  Implemented  for  Error  Parameters  in  SAIMUN 

is  the  true  angular  velocity  of  the  IMU  platform  with  respect  to  the  inertial  frame, 
expressed  in  the  IMU  frame  at  t^, 

TJMU  is  the  rotation  matrix  representing  the  rotation  from  the  inertial  to  the  IMU  frame, 
b()j.  is  the  startup  bias  of  the  accelerometers, 

ba  k  is  the  bias  of  the  accelerometer  at  tk,  which  changes  due  to  bias  instability, 

w0)k  is  the  thermo-mechanical  zero-mean  white  noise  present  in  the  accelerometers, 

Sa  is  the  scale  factor  error  matrix  of  the  accelerometers, 

Ma  is  the  axes  misalignment  matrix  of  the  accelerometers, 

Na  is  the  axes  nonorthogonality  matrix  of  the  accelerometers,  and 

aQ  k  is  the  quantization  effect  caused  by  analog-to-digital  conversion. 

and  similarly  for  the  gyroscopes.  For  simplicity,  let  alng  k  be  and  be  w*,.  The  error  sources 

are  applied  in  the  following  order: 

1.  Walking  bias  and  thermomechanical  noise  are  applied  first  because  they  affect  the  sensor 
(accelerometer  or  gyroscope)  regardless  of  how  the  sensor  is  mounted  with  respect  to  the 
defined  IMU  frame, 

2.  A  scale  factor  error  is  applied  next  to  account  for  errant  voltages,  circuitry,  etc.  in  converting 
the  sensor  output  to  a  value  that  can  be  quantized, 

3.  Axes  nonorthogonality  and  misalignment  errors  are  applied  next  to  account  for  the  mounting 
error  between  the  sensors  and  the  defined  IMU  frame, 

4.  Quantization  error  is  applied  last  to  emulate  the  Analog  to  Digital  Conversion  necessary  for 
quantizing  the  sensor  signal. 

The  mean  of  the  startup  bias  of  the  sensor  is  included  in  the  walking  bias  as  it  is  assumed  known 
from  sensor  testing. 

The  sources  of  error  used  by  this  model  affect  the  gyroscopes  identically  and  thus  their  presen¬ 
tation  is  omitted. 
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3.2.1  Startup  Bias.  The  startup  bias  affecting  the  accelerometers,  bay,  is  the  average  output  of 
the  accelerometers  while  they  are  not  undergoing  any  acceleration  [1],  For  practical  applications, 
the  startup  bias  is  quantified  by  accounting  for  and  removing  the  accelerometer  output  due  to 
gravity  while  the  IMU  remains  motionless.  Figure  (2)  shows  an  example  of  the  startup  bias  of  an 
accelerometer  undergoing  no  acceleration.  SAIMUN  assumes  that  the  statistics  of  the  startup  bias 


o 
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Figure  2.  Accelerometer  Output  Corrupted  by  Startup  Bias 

are  well  known  and  can  be  approximated  by  one  of  the  three  supported  error-parameter  pdfs.  The 
startup  bias  is  then  realized  from  the  error-parameter  pdf  during  virtual  IMU  initialization. 

3.2.2  Bias  Instability.  The  bias  of  an  accelerometer  tends  to  drift  slowly  over  time  and  is 
known  as  a  walking  bias  [1].  The  walking  bias  affecting  the  accelerometers  at  tk,  ba,fc>  is  the  change 
in  the  bias  after  the  IMU  was  powered  on.  This  effect  is  modeled  by  a  first-order  random  walk 
model,  which  causes  a  linear  increase  in  the  variance  of  the  accelerometer  bias  with  time.  The 
discrete  first-order  random  walk  model  is  given  by 

ba,k  =  ba<k_ i  +  wBl  k  ,  (2) 

where  wB1k  is  a  zero- mean  white  noise  process  of  standard  deviation  erBI.  The  initial  standard 
deviation  in  the  walking  bias  is  zero,  and  the  standard  deviation  at  a  later  time,  tBI,  is  specified 
as  tTgj.  Utilizing  the  zero  initial  condition  and  applying  Eq.  (2)  recursively  from  to  to  yields  the 
walking  bias  at  tk  as 


tk/At 

ba,k  —  ^  '  wBl,k  >  (3) 

i= 1 

where  At  is  the  timestep  of  the  IMU.  Noting  that  the  variance  of  the  summation  of  n  independent 
random  variables  denoted  by  Xj,  is  equivalent  to  the  sum  of  the  variances  of  n  independent  random 
variables,  i.e. 


this  yields 


Var 


V  Vai  V]  ; 
i= 1 


tk/At  t 

Var  [ba,k]  =  Var  =  i  ’ 

1=1 


(4) 
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when  applied  to  Eq.  (3).  Evaluating  Eq.  (4)  at  the  boundary  condition  gives 


=  a 


*  2 
BI  ' 


(5) 


The  standard  deviation  of  the  zero- mean  white  noise  process  driving  the  random  walk  model,  <xBI, 
can  now  be  found  by  manipulating  Eq.  (5)  as 


(6) 

Figure  (3)  shows  an  example  accelerometer  undergoing  a  random  walk  due  to  bias  instability. 
SAIMUN  inputs  the  boundary  conditions  <xBI  and  fBI,  the  IMU  timestep,  At,  and  the  desired  error- 


Figure  3.  Accelerometer  Output  Corrupted  by  Bias  Instability 

parameter  pdf-type  for  and  calculates  the  standard  deviation  of  the  underlying  noise  process 
driving  the  random  walk  model  according  to  Eqs.  (2)-(6).  Each  step  in  the  random  walk  model  is 
realized  independently  for  each  IMU. 

3.2.3  Velocity  Random  Walk.  When  thermo-mechanical  zero-mean  white  noise,  wa^,  present 
in  the  accelerometers  is  integrated,  a  random  walk  is  produced  in  the  velocity  state,  and  the 
standard  deviation  of  the  velocity  grows  proportionally  to  the  square-root  of  time  [1].  The  variance 
of  the  white  noise  present  in  the  accelerometers  is  described  by  the  Velocity  Random  Walk  (VRW) 
specification.  The  standard  deviation  of  the  zero-mean  white  noise  is  given  as  a  function  of  the 
VRW  specification  as 

_  VRW  spec 

^VRW  —  r 7 — 

\/A t 

Example  accelerometer  output  corrupted  by  zero-mean  white-noise  is  shown  in  Fig.  (4).  SAIMUN  in¬ 
puts  the  VRW  specification  and  IMU  timestep  to  calculate  the  standard  deviation  of  the  white  noise 
process  underlying  the  velocity  random  walk.  This  zero-mean  white  noise  process  is  realized  inde¬ 
pendently  at  each  tf;  for  the  virtual  IMU. 

3.2.4  Scale  Factor  Error.  Errors  in  analog-to-digital  converters  and  other  circuitry  can  induce 
scale  factor  errors  in  the  accelerometers,  denoted  by  Sa ■  For  instance,  errors  may  be  incurred  due 
to  reference  voltage  errors  caused  by  a  variety  of  sources,  including,  among  other  things,  tempera¬ 
ture  fluctuations  and  calibration  errors.  For  illustration  purposes,  consider  a  simple  example  of  a 
generic  scale  factor  error,  s,  and  assume  that  the  nominal  reference  voltage  for  the  analog-to-digital 


7 

Approved  for  public  release;  distribution  is  unlimited. 


Time  [s] 

Figure  4.  Accelerometer  Output  Corrupted  by  Thermo-Mechanical  Zero-Mean  White 
Noise 

converter  is  Fnom  =  5  V,  while  the  actual  reference  voltage  provided  to  the  analog-to-digital  con¬ 
verter  is  14ct  =  5.02  V.  The  scale  factor  error  of  the  signal  being  converted  by  the  analog-to-digital 
converter  is  then  given  by 

Fnnm  —  U,rt  5  V  —  5.02  V 

s  =  - —  = -  =  -0.003984  =  -0.3984% . 

Fact  5.02  V 

Denoting  the  scale  factor  errors  for  each  axis  of  the  accelerometer  triad  by  sx,  sy,  and  sz,  the  total 
accelerometer  scale  factor  error  may  be  expressed  as 

sx  0  0 

Sa  =  0  sy  0  . 

0  0  sz 

SAIMUN  inputs  an  error-parameter  pdf  for  sx,  sy,  sz  and  realizes  these  parameters  for  the  virtual 
IMU  from  this  pdf. 

3.2.5  Axes  Misalignment  Error.  Errors  in  the  manufacturing  process  and  mounting  of  an 
IMU  can  lead  to  axes  misalignment  errors,  denoted  by  Ma.  That  is,  the  actual  sensor  axes  are 
offset  from  the  IMU  case  frame  by  a  set  of  small  rotations.  This  is  readily  seen  by  considering  the 
term  I  +  Ma  in  Eq.  (la),  which  is  of  the  form  of  a  rotation  matrix  for  small  angles.  Therefore, 
the  actual  sensed  outputs  of  the  IMU  are  not  aligned  perfectly  with  the  modeled  case  frame  of  the 
IMU.  The  accelerometer  misalignment  errors  are  represented  within  Ma  by  rnx ,  my,  and  rriz ,  with 
each  term  corresponding  to  its  respective  sensor  axis.  The  misalignment  errors  are  expressed  in 
matrix  form  as 

0  mz  —iTiy 
Ma  =  —mz  0  mx 
my  —mx  0 

SAIMUN  inputs  an  error-parameter  pdf  for  rnx ,  my,  mz  and  realizes  these  parameters  for  the 
virtual  IMU  from  this  pdf. 

3.2.6  Axes  Nonorthogonality  Error.  Errors  in  the  manufacturing  process  of  an  IMU  can 
lead  to  axes  nonorthogonality  errors,  denoted  by  Na.  If  the  sensor  had  no  nonorthogonality  errors, 
each  axis  of  the  sensor  triad  would  be  at  perfect  right  angles  to  the  other  two  axes.  Therefore, 
the  nonorthogonality  errors  account  for  small  deviations  from  perfectly  orthogonal  axes  and  allow 
for  the  modeling  of  sensor  outputs  that  effectively  measure  small  overlaps  between  the  individual 
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axes.  The  accelerometer  axes  nonorthogonality  errors  are  represented  within  Na  by  nx,  ny,  and 
nz,  with  each  term  corresponding  to  its  respective  sensor  axis.  The  nonorthogonality  errors  can  be 
expressed  in  matrix  form  as 


Na 


0 

nz 

ny 

nz 

0 

rh 

ny 

nx 

0 

SAIMUN  inputs  an  error-parameter  pdf  for  nx,  ny ,  nz  and  realizes  these  parameters  for  the  virtual 
IMU  from  this  pdf. 

Example  output  of  an  accelerometer  triad  undergoing  1  m/s  acceleration  in  the  z  direction  with 
scale  factor,  axes  misalignment,  and  axes  nonorthogonality  errors  is  shown  in  Fig.  (5). 


s 
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Time  [s] 


Figure  5.  Three- Axis  Accelerometer  Output  Corrupted  by  Scale  Factor,  Misalignment, 
and  Nonorthogonality  Errors 


3.2.7  Quantization  Error.  Quantization  error  arises  in  IMUs  due  to  the  limitations  of  analog- 
to-digital  signal  conversion  and  is  modeled  such  that  the  acceleration  and  angular  velocity  mea¬ 
surements  are  affected  as 

(2n~ 1 

aQ ,k(ak)  =  round  (  - ak 

\  ®max 

where  n  is  the  bitrate  of  the  analog  to  digital  converter,  and  amax  is  the  maximum  value  able  to 
be  quantized  by  the  analog  to  digital  converter.  It  is  noted  that  the  quantization  error  is  applied 
after  all  of  the  other  errors  have  been  applied,  such  that  the  output  of  the  quantization  error  is 
the  measured  acceleration  of  Eq.  (la).  Values  of  the  bitrate  and  maximum  value  are  both  input 
into  SAIMUN.  The  quantization  error  in  angular  velocity  is  calculated  similarly.  Figure  (6)  shows 
a  quantized  signal  from  an  IMU  with  thermo- mechanical  noise  present. 


®max 
2n~1  ’ 
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Figure  6.  Accelerometer  Output  Corrupted  by  Quantization  Error 


3.3  IMU  Inversion 

In  order  to  perform  dead-reckoning  navigation,  the  sensor  errors  must  be  removed  from  the  non- 
gravitational  acceleration  and  angular  velocity.  An  expression  for  how  these  quantities  are  corrupted 
by  the  sensor  are  given  by  Eqs.  (1).  Solving  this  equation  for  the  uncorrupted  measurements  requires 
an  “inversion.”  The  measured  non-gravitational  acceleration  and  angular  velocity  are  modeled  using 
Eqs.  (1).  For  some  vector  v  =  [vx  vy  vz]T ,  define  the  matrices  [u\],  [ux],  and  [u*]  to  be 


Vx 

0 

0 

0 

Vz 

- Vy 

0 

Vz 

Vy  ' 

v\]  = 

0 

Vy 

0 

,  [vx]  = 

~Vz 

0 

Vx 

,  and  [v*]  = 

vz 

0 

Vx 

0 

0 

vz  _ 

.  vy 

Vx 

0  j 

.  vv 

Vx 

0 

After  omitting  the  effects  due  to  quantization,  the  errors  in  the  acceleration  and  angular  velocity 
can  be  expressed  as 

am,k  =  {I  +  Na  +  -/Via)  (7  +  Sa )  (ttfc  +  ba +  Wa^)  ,  (7) 

where  a &  =  T’/MUai  is  used  for  compactness.  Equation  (7)  may  be  solved  for  a ^  in  terms  of  the 
measured  acceleration  and  the  error  sources  to  yield 

d’k  =  (I  T  So)  ( I  +  Na  +  Ma)  CLm,k  —  ba,k  VJ a,k  ■ 

Noting  that  ( 7  +  Ma  +  Na)  (7  +  Sa )  ~  7  +  Aa,  applying  the  matrix  inversion  lemma,  and  sim¬ 
plifying  the  resulting  expression,  it  follows  that  the  true  non-gravitational  acceleration  written  in 
terms  of  the  measured  non-gravitational  acceleration  as 

T  [®m,fcX]ma  ^a,fc  V) a,k  ;  (8) 

from  which  one  may  obtain  an  estimate  of  the  true  non-gravitational  acceleration  as  a =  Eja^}, 
which  gives 

dfc  T  [dr n,k*\^a  ^a,fc  ^a,k  •  (9) 

If  all  of  the  error  sources  are  zero-mean,  it  follows  from  Eq.  (9)  that 

U'k  —  0'm,k  ■ 

Parallel  results  hold  for  expressing  the  true  angular  velocity  in  terms  of  the  measured  angular 
velocity.  That  is,  following  the  same  process  used  in  arriving  at  Eqs.  (8)  and  (9),  it  can  be  shown 
that  the  true  and  estimated  angular  velocities  can  be  expressed  as 

^ k  =  m,k  T  \uJrn1k*\^,g  ~  —  ^ g,k  IV g,k 

&k  =  W. m,k  T  ~  [(Vm,k*]'R/g  —  —  bg,k  IV g,k 
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where  all  of  the  error  sources  are  now  for  the  gyro  instead  of  the  accelerometer.  As  with  the 
accelerometer,  if  all  of  the  error  sources  are  zero-mean,  it  follows  that 

^ k  =  W m,k  ■ 

Define  the  non-gravitational  acceleration  and  angular  velocity  deviations  to  be  the  difference  in 
the  truth  from  the  estimate  as 


d  CLk  —  ft  ft 

5u:k  =  uk  -  uk  . 

Taking  the  difference  between  true  and  estimated  expressions  yields 

hftfc  =  X  \o,mk^\8na  [ft'jn,fcX']<^^a  3ba,k  &U) ak 

=  —  —  ^ g,k  ~  $^g,k  j 

where  the  deviations  in  the  error  parameters  are  defined  similarly  to  acceleration  and  angular 
velocity. 

3.4  IMU  Mechanizations 

Traditionally,  IMU-based  navigation  is  mechanized  in  one  of  two  ways:  gimbaled/gyro-stabilized  or 
strapdown.  In  the  gimbaled  mechanization,  the  triad  of  accelerometers  is  placed  on  a  “stable  table” 
at  the  center  of  set  of  interconnected  gimbals.  This  enables  the  stable  table,  through  the  use  of 
feedback  tracking  using  the  gyros,  to  track  a  predetermined  attitude  such  that  the  accelerometers 
measure  a  certain,  desirable  resultant  linear  motion.  Typically,  the  stable  table  is  configured  to 
track  a  north-east-down  coordinate  system,  in  which  case  the  accelerometers  measure  the  linear 
motion  in  these  directions,  or  it  is  configured  to  track  the  inertial  coordinate  system  (i.e.  to  remain 
fixed  in  space  as  the  vehicle  maneuvers). 

The  most  basic  of  gimbal-based  mechanizations  employs  three  interconnected  girnbal  rings, 
such  as  is  found  on  the  Minuteman  III  intercontinental  ballistic  missile  (see  Fig.  (7(a))).  While 
this  configuration  provides  the  ability  to  actively  control  the  orientation  of  the  stable  table  through 
motors  mounted  on  the  gimbals,  the  use  of  three  gimbals  sometimes  leads  to  girnbal  lock  and 
prevents  the  stable  table  from  tracking  all  possible  attitudes  that  may  be  required. 

More  advanced  configurations  of  the  stable  table  mounting  design  can  be  used  to  replace  the 
original  gimbaled  system.  This  is  accomplished  by  using  four  or  more  gimbals  or  by  using  a  fluid 
stabilization  system,  which  has  the  stable  table  floating  in  the  center  of  a  sphere  filled  with  pressur¬ 
ized  fluid.  Whereas  the  attitude  of  the  stable  table  is  controlled  using  torquers  and  manipulation 
of  the  gimbals  in  the  original  design,  the  attitude  of  the  stable  table  in  the  fluidic  design  is  con¬ 
trolled  via  nozzle  actuators  which  avoids  any  possible  complications  due  to  girnbal  lock.  One  such 
implementation  of  the  fluidic  system  is  the  Advanced  Inertial  Reference  Sphere  (AIRS)  that  was 
developed  for  the  LGM-118A  Peacekeeper  (see  Fig.  (7(b))). 

In  the  strapdown  mechanization,  the  triads  of  accelerometers  and  gyros  are  fixed  (strapped 
down)  to  the  body  of  the  vehicle,  and  there  is  no  stable  table.  This  means  that  where  the  stable 
table  was  used  to  mechanically  track  the  inertial  frame  in  order  to  track  a  predetermined  attitude 
and  sense  accelerations  in  this  rotating  frame,  the  gyro  outputs  are  now  used  to  numerically  track 
the  orientation  of  the  IMU  (called  the  case  frame)  with  respect  to  the  inertial  frame.  To  avoid  any 
singularities  present  in  three-parameter  attitude  representations  such  as  Euler  angles,  the  attitude 
quaternion  is  used  and  can  be  propagated  forward  in  time  using  the  sensed  angular  velocity  from 
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(a)  Minuteman  III  IMU  (b)  Peacekeeper  IMU 

Figure  7.  Examples  of  Stable  Table  IMU  Mechanizations  [2] 


the  IMU.  Additionally,  the  accelerations  are  now  sensed  in  the  case  frame,  and  these  must  be 
transformed  using  the  case-to-inertial  attitude  in  order  to  predict  the  translational  position  and 
velocity  of  the  vehicle. 

Two  inertial  navigation  system  mechanizations  are  presented  here:  a  strapdown  system  and 
a  passively  commanded  inertially  stable  space  stabilized  gimbaled  system.  For  each  of  the  afore¬ 
mentioned  mechanizations,  the  governing  equations  for  the  position,  velocity,  and  attitude  of  the 
IMU  are  presented,  and  a  set  of  discretized  equations  are  developed.  The  discretized  equations 
governing  the  position,  velocity,  and  attitude  of  the  IMU  then  are  representative  of  the  software 
and  hardware  that  would  be  used  in  each  inertial  navigation  mechanization. 

In  the  following  sections,  the  following  remains  consistent: 


v 

q 

a 

Ul! 

U) 

rpb 

-*■  n 


position 

velocity 

attitude  expressed  as  a  quaternion 

acceleration 

angular  velocity 

angular  velocity  expressed  as  a  true  quaternion 
rotation  matrix  from  frame  a  to  frame  b 


IMU 

CM 

P 

b 

s 

i 

f 

9 

ng 

m 


inertial  measurement  unit 
center-of-mass  of  the  vehicle 
platform  frame  of  the  IMU 
body  frame  of  the  vehicle 
sensor  case  frame 
Earth-Centered  Inertial  frame 
Earth-Centered  Earth-Fixed 
due  to  gravity 
not  due  to  gravity 
measured 


and  a  subscript  given  in  the  form  a/b  denotes  a  with  respect  to  b.  For  brevity,  only  the  governing 
equations  and  the  resulting  discretized  update  equations  are  given  in  the  following  sections. 
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3.4.1  Strapdown  Mechanization.  The  translational  and  rotational  states  of  a  strapdown  IMU 
aboard  a  vehicle  are  governed  by  [3,  4] 

^imu(^)  =  ®imuW 

«i'mu  (*)  =  al(rl  mu(0  +  Tb(t)r  cm/imu)  +  Tb(t)ang(t) 

h\it)  =  \ubb/iit)  ®  Qiit) , 

For  notational  ease,  let 


riMV(t) 

— >• 

r{t)  , 

^IMU  it) 

+  v(t)  , 

<(•) 

— )• 

9ih 

mt) 

-» 

TT(t)  , 

abng(t)  ~ 

>  ait)  , 

Qi.it) 

— > 

Qit) 

* b/iit ) 

— >• 

“it)  , 

rp 

cm/imu 

d  , 

and  rtMV(t)  +  T£it)rbCM/lMV 

-4 

s{t) 

With  these  substitutions,  the  equations  of  motion  may  be  written  succinctly  as 

f{t)  =  v(t) 

v(t)  =  g(s(t))  +Tr(t)a(t) 
hit)  =  >{t)®q(t). 

The  discrete-time  dead-reckoning  equations  for  the  position,  velocity,  and  attitude  of  a  vehicle 
with  a  strapdown  IMU  are  given  by: 

rk  =  rfc_  i  +  •Ufc-iA  tk  +  ^gk-iAtl  +  ^Tk_iakAtk  (10a) 

—  g  iGk-iTk_i[dx]  +  Tk_\  [ofcx])  ukAt\ 

vk  =  ufc_i  +  gk-iAtk  +  T^_xakAtk  -  *  (G^T^dx]  +  Tk_x  [akx])  u>kA t\  (10b) 

qk  =  q(ukAtk)  <8)  qfc_ i ,  (10c) 

where  gk-\  and  Gk-\  are  the  two-body  gravitational  acceleration  experienced  at  the  CM  (not  the 
IMU)  and  the  associated  Jacobian  evaluated  at  the  estimated  position  of  the  CM  at  time  tk_\. 

3.4.2  Space  Stabilized  Mechanization.  The  space  stabilized  mechanization  tracks  an  Earth- 
centered,  inertially  non-rotating  coordinate  frame.  This  reference  frame  is  realized  by  a  gyro- 
stabilized  platform  with  triaxial  control.  Except  for  torques  applied  as  compensation  to  account 
for  anisoelastic  effects,  thermal  and  acceleration  sensitivity,  etc.,  the  gyroscopes  are  uncommanded 
and  the  system  is  “free-floating.” 

Noting  that  the  specific  force  is  actually  measured  in  the  IMU  platform  frame,  the  transla¬ 
tional  equations  governing  the  dynamical  evolution  of  the  position  and  attitude  of  the  IMU  can  be 
expressed  as  a  set  of  first-order  differential  equations  given  by 

^IMU  (t  )  —  ^IMU  it) 

<U  it)  =  9lir\MV{t ))  +  Tp(t)fp(t) 

hi  it)  =  \upp/i{t)®  qp(t) 
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For  notational  simplicity,  let  rl1MV(t)  ->•  r,  vl1MV(t)  v,  upp/i{t)  u(t),  qf{t)  q(t ),  T*(t)  = 

Tt (t),  and  f  p(t)  =  f(t),  such  that  the  governing  equations  become 

r(i)  =  v(t) 

v{t)  =  g{r(t))  +  TT(t)f{t) 
q{t)  =  ^u)(t)®q(t). 

Discretization  and  integration  yield  the  following  discrete-time  propagation  equations: 


rk  =  r*fc_i  +  4  +  ^  {Tk_ifk  +  gk- 1)  A tk 

(11a) 

vk  =  +  (T^LjOfc  +  gk- 1)  At.k 

(lib) 

qk  =  q(u)kAtk)  <8)  gfc_i . 

(11c) 

3.5  Monte  Carlo  Analysis 

In  general,  either  of  the  mechanizations  of  the  IMU  can  be  represented  by  a  nonlinear  dynamical 
system.  If  the  expressions  for  position,  velocity,  and  attitude  are  collected  into  a  state  vector,  x , 
the  general  nonlinear  dynamical  system  representing  the  forward  evolution  of  the  state  for  either 
of  the  IMU  systems  may  be  expressed  as 

f  (®fc—  1 )  ^ki  U)  j  (12) 

where  ak  and  ujk  represent  the  true  vehicle  non-gravitational  acceleration  and  angular  velocity, 
respectively.  For  the  mechanizations  presented  in  Section  3.4,  /(•)  represents  one  of  Eqs.  (10)  or 
Eqs.  (11).  Functionally,  the  relationship  between  the  true  and  measured  IMU  outputs  is  given  by 

am,k  =  a(<f)k,tk )  and  =  oj(i/)k,tk) ,  (13) 

where  4>k  and  ipk  represent  the  IMU  model  parameters,  which  also  have  governing  evolutionary 
equations  as 


4>k  =  9(<t>k,tk )  and  tpk  =  h(ij>k,tk) .  (14) 

In  the  preceding  discussions,  many  possibilities  for  the  forms  of  g ,  h ,  <f>,  and  -0  have  been  given.  For 
example,  for  the  simple  additive  model  where  the  true  IMU  outputs  are  corrupted  by  a  constant 
bias  b  and  a  white- noise  sequence  w ,  Eqs.  (13)  become 

am,k  =  ak  +  ba  +  wa)k  and  <^m,k  =  <^k  +  bg  +  wg)k  , 

and  the  parameters  of  the  IMU  error  model  are 


0/c  — 

ba 

and 

= 

*>9 

_  wa,k  _ 

_  Wg,k  _ 

with  the  governing  equations  cf)k  =  (f)k-\  and  tpk  =  xjjk- 1  (note  that  the  white  noise  sequence  does 
not  require  a  formal  time- wise  evolution  owing  to  the  properties  of  white  noise). 

Given  the  description  of  the  IMU-based  dynamical  system  via  Eq.  (12),  the  dynamical  descrip¬ 
tion  of  the  IMU  model  parameters  via  Eq.  (14),  and  a  statistical  description  of  the  IMU  model 
parameters  that  appear  in  Eq.  (13),  the  objective  is  to  determine  a  statistical  description  of  the 
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position,  velocity,  and  attitude  of  the  IMU  at  some  later  time.  While  there  is  no  tractable  exact 
solution  for  predicting  these  statistics,  in  general,  several  approaches  exist  for  approximating  the 
statistics.  One  method  is  the  Monte  Carlo  process.  The  Monte  Carlo  process  effectively  performs 
a  large  number  of  simulations  of  the  dynamical  system,  with  each  simulation  drawing  new  sam¬ 
ples  from  every  statistical  distribution  involved  in  the  system.  As  such,  the  accuracy  of  a  Monte 
Carlo  process  is  limited  by  the  number  of  simulations  considered,  which  is,  in  turn,  limited  by 
the  computational  time  allotted.  On  the  other  hand,  the  Monte  Carlo  process  provides  a  very 
flexible  environment  that  is  capable  of  handling  a  wide  variety  of  dynamical  systems  and  statistical 
distributions. 

To  illustrate  the  Monte  Carlo  process,  consider  the  two  dimensional  transformation  from  polar 
to  Cartesian  coordinates,  i.e. 


x 

V 


=  r 


cos  9 
sin0 


Given  a  distribution  on  w  =  [r  8\,  it  is  desired  to  determine  distribution  characteristics  on 
z  =  [x  y] .  First,  assume  that  r  and  9  are  independent  and  Gaussian  with  mean  and  covariance 

(0.02  [m])2  0 

0  (30  [deg])2  _ 

Therefore,  the  polar  to  Cartesian  conversion  may  be  written  as 


mw  = 


1  [m] 
60  [deg 


and 


Pw  = 


z  =  f{w) , 

where  w  is  the  input  and  z  is  the  associated  output.  The  Monte  Carlo  simulation  process  is 
started  by  drawing  a  number  of  samples  from  the  input  vector,  where  the  sampling  is  performed 
with  respect  to  the  probability  density  function  which  governs  the  characteristics  of  the  input. 
Once  the  samples  are  drawn,  the  nonlinear  mapping  function,  /,  is  applied  to  each  input  sample 
to  generate  a  set  of  output  samples. 

For  the  coordinate  conversion  problem,  the  Monte  Carlo  simulation  is  initialized  by  drawing 
1  X  105  realizations  of  the  random  vector  w  from  the  distribution 

p(w;mw,Pw )  =  |27rPwp1/2exp  |-^  ( w  -  mw )T  P”1  ( w  -  mw) 

The  set  of  samples  drawn  is  shown  in  Fig.  (8(a)).  The  nonlinear  transformation  is  applied  to  each 
sample  in  the  input  space  and  the  resulting  set  of  Cartesian  coordinates  is  shown  in  Fig.  (8(b)).  It 
is  clear  that  the  output  space  is  not  capable  of  being  described  by  a  Gaussian  distribution;  however, 
the  Monte  Carlo  simulation  technique  can  still  be  used  to  compute  an  accurate  approximation  of 
the  mean  and  covariance  via 

1  N 

and  Pz  = —^(zi  -  mz)(zi  -  mz)T  ,  (15) 

1=1 


1 


N 


mz  =  N^Zl 

i= 1 


where 


Zi  =  f(Wi) 

and  Wj  is  the  ith  Monte  Carlo  sample.  The  contours  representing  the  Gaussian  distribution  com¬ 
puted  from  the  preceding  mean  and  covariance  equations  are  shown  on  top  of  the  set  of  output-space 
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r  [m]  x  [m] 

(a)  Input  Space  Monte  Carlo  Samples  (b)  Mapped  Samples  in  the  Output  Space 

Figure  8.  Input  and  Output  Samples  from  the  Monte  Carlo  Simulation 


Monte  Carlo  samples  in  Fig.  (9(a)).  Moreover,  a  linearization  approach  can  be  used  to  map  the 
input  mean  and  covariance  into  an  approximate  output  mean  and  covariance.  The  comparison 
between  the  linear  mean  and  covariance  and  the  Monte  Carlo  mean  and  covariance  is  shown  in 
Fig.  (9(b)),  wherein  it  is  seen  that  the  two  differ  a  great  deal.  The  only  characteristic  that  is 
mostly  captured  by  the  linearization  process  is  the  orientation  of  the  covariance  matrix.  The  cen¬ 
tering  (i.e.  the  mean  of  the  distribution)  is  not  captured  accurately,  and  the  size  of  the  mapped 
covariance  matrix  is  not  captured  accurately. 


i 

a 

0.5 


(a)  Monte  Carlo  Mean  and  Covariance  (b)  Comparison  between  Linearization  and  Monte  Carlo 

Figure  9.  Distribution  Contours  Obtained  by  the  Mean  and  Covariance  from  the 
Monte  Carlo  Simulation  and  Linearization 

To  further  demonstrate  the  power  of  Monte  Carlo  simulation,  consider  a  modification  of  the 
preceding  problem.  In  this  case,  it  is  still  assumed  that  the  radial  and  angular  input  coordinates 
are  independent,  but  now  the  radial  direction  is  taken  to  be  Gaussian  distributed  and  the  angular 
direction  is  taken  to  be  uniformly  distributed.  The  mean  and  covariance  are  taken  to  be  the  same 
as  before  so  that  the  input  space  has  the  same  first  and  second  central  moments  as  in  the  preceding 


x  [ml  x  [ml 


16 

Approved  for  public  release;  distribution  is  unlimited. 


example.  The  probability  density  function  for  the  input  is  now  given  by 


p(w;mr,Pr,mo,Pg) 


1 2irPr  |  ^2  exp 

2V6 TrPrPo  CXP 


(wr  —  mr )2  1  1 

2Pr  J  '  2VWe 
(wr  —  mr )2  | 

2Pr  j  ’ 


where  me  —  \/3Pg  <  6  <  mg  +  \/3 Pg,  such  that  [mg  —  \/3 Pg  ,  mg  +  \/3Pg\  is  the  support  set  of  the 
uniform  angular  distribution.  From  the  original  specifications  for  the  transformation  problem,  the 
means  and  covariances  required  to  evaluate  the  input-space  distribution  are 


mr  =  1  [m]  ,  mg  =  60  [deg]  ,  Pr  =  (0.02  [m])2  ,  and  Pg  =  (30  [deg])2  . 

As  before,  the  Monte  Carlo  simulation  is  initialized  by  drawing  1  x  105  realizations  of  the  random 
vector  w  from  the  distribution  p(w;mr,  Pr,mg,  Pg),  and  each  sample  is  mapped  into  the  output 
space  via  the  nonlinear  transformation.  These  samples  are  shown  in  Figs.  (10(a))  and  (10(b)), 
respectively.  The  mean  and  covariance  is  again  computed  from  the  Monte  Carlo  samples  and 
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(a)  Input  Space  Monte  Carlo  Samples 
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(b)  Mapped  Samples  in  the  Output  Space 


Figure  10.  Input  and  Output  Samples  from  the  Monte  Carlo  Simulation  (modified 
problem) 

shown  along  with  the  output-space  samples  in  Fig.  (11(a)).  Finally,  a  linearization  approach  is 
used  to  map  the  input-space  mean  and  covariance  into  an  output-space  mean  and  covariance. 
The  comparison  between  the  Monte  Carlo  mean  and  covariance  and  the  linearization  mean  and 
covariance  is  shown  in  Fig.  (11(b)). 

While  it  is  clear  that  the  linearization  process  once  again  fails  to  capture  the  true  mean  and 
covariance,  it  is  noted  that  the  linearized  mean  and  covariance  depicted  in  Fig.  (11(b))  are  identical 
to  the  ones  depicted  in  Fig.  (9(b)).  At  the  same  time,  however,  the  Monte  Carlo  means  and 
covariances  differ  between  these  two  simulations.  This  occurs  because  the  Monte  Carlo  simulation 
has  and  uses  knowledge  of  the  entire  distribution  on  both  the  input  and  output  spaces  whereas  a 
linearization  procedure  works  with  only  the  first  and  second  central  moments  of  the  distribution. 
Therefore  if  two  input  distributions  have  the  same  mean  and  covariance,  and  are  subjected  to  the 
same  nonlinear  transformation,  the  linearized  solutions  will  always  have  the  same  output  mean  and 
covariance. 
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1.2 


(a)  Monte  Carlo  mean  and  covariance 


(b)  Comparison  between  linearization  and  Monte  Carlo 


Figure  11.  Distribution  Contours  Obtained  by  the  Mean  and  Covariance  from  the 
Monte  Carlo  Simulation  and  Linearization  (modified  problem) 


It  is  also  worth  noting  that  the  Monte  Carlo  simulation  enables  analysis  extending  beyond  the 
computation  of  the  mean  and  covariance.  For  instance,  higher-order  statistical  moments  may  be 
computed  in  a  similar  fashion  to  the  mean  and  covariance,  best  case  and  worst  case  performance 
can  be  analyzed  through  analysis  of  the  dispersions  away  from  a  predefined  truth,  or  the  radial 
distance  encapsulating  a  given  percentage  of  the  sample  points  can  be  computed. 

The  problem  of  analyzing  the  performance  of  IMU-based  navigation  systems  is  characterized 
by  complex,  nonlinear  functions  where  the  input  space  is  described  by  the  IMU  error  sources, 
which  may  have  fairly  arbitrary  statistical  descriptions.  Monte  Carlo  analysis  naturally  handles 
the  characteristics  of  this  problem  and  enables  a  more  rich  data  analysis  suite  than  is  found  in 
typical  linearization-based  analysis. 

3.6  External  Aiding 

In  this  section,  models  for  three  different  external  sensors  available  within  SAIMUN  are  described. 
The  sensors  presented  are 

1.  range/range-rate, 

2.  stellar  line-of-sight,  and 

3.  bearing  angles. 

For  each  sensor  considered,  the  measurement  model  is  described  for  each  of  the  previously  described 
IMU  mechanizations.  Additionally,  to  facilitate  data  processing,  the  measurement  Jacobian  for  each 
configuration  is  described. 

3.6.1  Range/Range- Rate  Measurements.  The  range  f>,  by  definition  is  the  difference  mag¬ 
nitude  between  two  position  vectors.  For  this  case,  the  two  position  vectors  will  be  the  vehicle’s 
position  in  a  given  frame  rs,  and  the  position  vector  of  a  transmission  receiving  ground  station  rg 
in  the  same  frame.  rg  is  assumed  to  be  known. 
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Strapdown  The  range  measurement,  zp,  is  modeled  as  the  range  between  an  onboard  sensor  and 
a  known  reference  position  that  is  corrupted  by  additive,  zero-mean,  white  measurement  noise;  the 
corresponding  model  is 

ZP  =  +  VP  > 

where  vp  is  the  measurement  noise  and  rs/g  is  the  true  position  of  the  range  sensor  with  respect  to 
the  known  ground  station,  which  is  given  by 

<h  =  <,v  +  Tlrb„lMV-T}  r>. 

A  reference  measurement  may  be  determined  by  evaluating  the  measurement  model  at  a  reference 
state,  leading  to  zp.  The  deviation  between  the  true  and  reference  measurements  is  then  considered 
as 


5zp  =  zp  —  z 


p  • 


This  deviation  is  then  linearized  using  a  first-order  Taylor  series  expansion,  and  the  result  may  be 
expressed  as 


$zp  —  HPir6r1MV  +  HpfibQ  +  vp  , 


where  <5rjMU  is  the  error  between  the  true  and  reference  IMU  positions  and  56  is  the  attitude  error, 
which  is  defined  to  be  twice  the  vector  part  of  the  quaternion  error.  This  illustrates  that  the 
range  sensor  measurement  model  has  sensitivity  associated  with  the  position  of  the  IMU  and  the 
attitude.  This  means  that  direct  updates  can  be  obtained  for  these  states  when  considering  range 
measurements  from  a  known  ground  station. 

Similar  to  the  range  measurement,  the  range-rate  measurement  produces  a  measurement  of  the 
range-rate  between  an  onboard  sensor  and  a  known  ground  station.  This  is  represented  by  zp, 
where 


zP  = 


(ris/g)T«/g) 


+  Vp, 


where  vp  is  zero  mean,  white,  range-rate  measurement  noise,  rs/g  is  the  position  between  the  sensor 
and  the  ground  station,  and  the  velocity  of  the  range-rate  sensor  with  respect  to  the  ground  station 
are 


v: 


Js/g  —  UIMU  T  ±b  VJJb/iA' \'  s/imv 

As  with  the  range  measurement,  a  reference  measurement  may  be  determined  by  evaluating  the 
measurement  model  at  a  reference  state,  leading  to  zp.  The  deviation  between  the  true  and  reference 
measurements  is  then  considered  as 


+  Ti[ubb/ix]rbs 


J 


f 


~  THU}f/iX\rg 


OZp  =  Zp  -  Zp. 

This  deviation  is  then  linearized  using  a  first-order  Taylor  series  expansion,  and  the  result  may  be 
expressed  as 

5zp  =  Hpr5v1MV  +  H P'V5vimv  +  HpM56  +  Hp.bg5bg  5mg  +  Hp^ng5ng  +  HP:Sg5sg  +  vp , 

where  is  the  error  between  the  true  and  reference  IMU  velocities,  5bg  is  the  gyro  bias  error, 

5mg  is  the  gyro  misalignment  error,  5ng  is  the  gyro  nonorthogonality  error,  and  5sg  is  the  gyro 
scale  factor  error.  Thus,  the  range-rate  measurement  has  sensitivity  associated  with  the  position 
of  the  IMU,  the  velocity  of  the  IMU,  the  attitude  of  the  vehicle,  and  the  gyro  parameters. 
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Space  Stabilized  The  range  measurement  for  the  space  stabilized  case  is  nearly  the  same  as  for 
the  strapdown  case;  therefore,  the  model  can  be  written  as 


^  =  \ArWT(rW+^’ 

where  vp  is  zero-mean,  white,  measurement  noise.  The  term  rs/g  still  represents  the  position  of  the 
range  sensor  with  respect  to  the  known  ground  station,  but  is  slightly  changed  from  the  strapdown 
case,  such  that  is  it  now  given  by 


s/g 


=  +  T'X rJ/IU„  -  Tirl 


1  9 


A  reference  measurement  may  be  determined  by  evaluating  the  measurement  model  at  a  reference 
state,  leading  to  zp.  The  deviation  between  the  true  and  reference  measurements  is  then  considered 


as 


6zp  =  zp-  zp. 


This  deviation  is  then  linearized  using  a  first-order  Taylor  series  expansion,  and  the  result  may  be 
expressed  as 


5zp  —  H p  r5vlMlj  T  Hp  o5d  +  vp  . 


It  is  important  to  note  that  the  attitude  error,  56,  which  is  still  defined  to  be  twice  the  vector  part 
of  the  quaternion  error,  now  defines  the  attitude  error  of  the  platform,  as  opposed  to  of  the  body. 

As  with  the  space  stabilized  range  measurement,  the  range-rate  measurement  model  is  largely 
unchanged  by  the  IMU  mechanization;  therefore,  the  model  may  be  written  as 


zP  = 


(rWT(*4) 


+  Vp, 


where  vp  is  zero  mean,  white,  range-rate  measurement  noise,  and  the  position  and  velocity  relative 
to  the  ground  station  are 


-r*  — 

s/g  in 


I  rp  i  rp  'P  b  _  rpi  f 

'J-p-Lb's/ imu  1  frg 


</a  =  +  TvT>  b/v  )  ^s/imu  f  H  ^ 


LP^p/i 


p  b  b/p' 


s  I  imu 


J 


The  deviation  between  the  true  and  reference  measurements  is  then  considered  as 


5zp  =  zp-  zp, 

where  zp  is  the  reference  measurement.  This  deviation  is  then  linearized  using  a  first-order  Taylor 
series  expansion,  and  the  result  may  be  expressed  as 

5zp  Hp)r5rlMV  +  H p'V5vlMXj  +  HpjpO  +  Hp^5bg  Hp,mg  5mg  +  Hp,ng5ng  +  HP:Sg8sg  +  vp  . 

As  with  the  range  measurement  for  the  space  stabilized  IMU  mechanization,  the  attitude  error, 
56,  represents  the  platform  attitude  error.  It  should  also  be  noted  that  the  exact  terms  appearing 
in  the  preceding  first-order  expansion  are,  in  general,  different  than  those  found  for  the  strapdown 
range-rate  measurement. 
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Data  Processing  When  considering  both  range  and  range-rate  measurements,  the  combined 
measurement  and  reference  measurement  are 


zp,p  ~ 


and 


zp,p  = 


In  order  to  process  these  measurements,  it  is  required  to  determine  the  Jacobian  of  the  combined 
measurement.  The  Jacobian  for  the  combined  measurements  is  simply  given  by  the  concatenation 
of  the  individual  Jacobians,  or 

Hp,r  03x3  Hpfi  03x3  O3X3  O3X3  O3X3  O3XI2 

HpQ  Hpfj  Hppn  Hp^n  Hpts  O3X9  J 


3.6.2  Stellar  Line-of-Sight  Measurements.  Star  tracking  measurement  systems  utilize  a 
camera  and  the  knowledge  provided  by  star  field  databases  to  assist  in  attitude  determination.  For 
simplicity,  the  measurements  considered  here  are  unit  vectors  to  each  star  in  the  star  camera  field 
of  view.  However,  this  solution  suffers  from  the  fact  that  the  covariance  matrix  of  a  unit  vector,  as 
well  as  the  addition  used  in  the  additive  update,  are  not  well  defined.  Two  alternative  approaches, 
focal  plane  angles  (FPA)  and  quaternion  estimation  (QUEST),  avoid  such  issues,  and  should  be 
considered  for  future  investigation. 


Strapdown  The  model  of  the  star  camera  measurement  is  given  as 


%SC 


_  rpscrpbrpl  sr 

~  1b  ±i  Jsrus/sc 


+  V 


sc 

sc 


where  Tbc  is  the  rotation  matrix  relating  the  body  and  star  camera  reference  frames,  T]J  is  the 
rotation  matrix  relating  the  inertial  and  body  frames,  T*r  is  the  rotation  matrix  relating  the  stellar 
reference  and  inertial  frames,  is  zero-mean,  white,  measurement  noise,  and  usJjsc  is  the  unit 
vector  representation  of  the  line-of-sight  to  a  reference  star,  which  can  be  expressed  in  terms  of 
right  ascension  ( asr )  and  declination  (e)sr)  as 


u 


s/ SC 


cos  5sr  cos  asr 
cos  5sr  sin  asr 
sin  5sr 


Similarly,  the  reference  star  camera  measurement  can  be  written  as 


~SC  _  rjiSCrpbrrii  „.ST 

zsc  srus/sc  * 


The  deviation  of  the  star  camera  measurement,  which  is  the  difference  between  the  true  and  refer¬ 
ence  measurements,  can  be  expressed,  after  a  first-order  Taylor  series  expansion,  as 

5z?c  =  Tbsc  [TfT>f/scx]  SO  +  vfc  , 

where  56  is  the  attitude  error,  and  it  is  noted  that  this  represents  the  attitude  error  for  the  body- 
fixed  reference  frame. 
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Space  Stabilized  The  space  stabilized  mechanization  model  is  nearly  identical  to  the  strapdown 
mechanization  model,  except  that  an  extra  rotation  must  be  introduced  to  account  for  the  ap¬ 
pearance  of  the  platform  frame.  Thus,  the  measurement  generated  by  the  star  camera  is  modeled 
as 


z 


SC 

SC 


rjiSCrrtbrriPrjii  Sr 

-Lb  -Lp-Li  -LsrUjs/sc 


+  v 


sc 

SC  5 


and  the  reference  measurement  is  found  via 


~SC  _  rr\SCrpbrpPrr\i  -.ST 

Zsc  —  ±b  1p±i  1srUs/sc • 

As  with  the  strapdown  case,  the  difference  between  the  measurement  and  reference  measurement  is 
defined  to  be  a  measurement  deviation,  5z^,  and  this  deviation  is  expanded  in  a  first-order  Taylor 
series  to  yield 


=  ncTbp  [TfTirus;/scy]  56  +  vfc  , 

where  56  is  the  attitude  error,  and  it  is  noted  that  this  represents  the  attitude  error  for  the 
platform-fixed  reference  frame. 


Data  Processing  In  order  to  process  stellar  line-of-sight  measurements,  it  is  required  to  deter¬ 
mine  the  Jacobian  of  the  measurement.  The  Jacobian  for  the  line-of-sight  measurements  takes  on 
the  same  form,  independent  of  the  IMU  mechanization,  such  that  it  can  be  written  as 

Hsc  =  [03x6  Hg  03x24]  , 


where 


He  =  Tbsc[fibHsrus;/scx\ 
for  the  strapdown  IMU  mechanization  and 

He  =  ncTbp\ffTirus;/s^] 

for  the  space  stabilized  IMU  mechanization. 


3.6.3  Bearing  Angles  Measurements.  A  coordinate  system  is  defined  such  that  the  y  axis 
points  along  the  “nose”  of  the  vehicle,  the  z  axis  points  “up,”  and  the  x  axis,  which  completes  the 
triad,  points  “right.”.  The  bearing  angles  measurements  are  measurements  of  the  line-of-sight  to  a 
destination  target,  where  the  angles  are  taken  to  be  the  elevation  and  azimuth  angles.  The  relative 
position  of  the  destination  target  with  respect  to  the  bearing  angles  sensor  is  denoted  by  pt/s  (or 
simply  p ,  for  ease  of  exposition),  where  the  destination  target  is  assumed  to  be  known  in  the  fixed 
frame.  The  elevation  angle  is  measured  on  the  interval  [— ir/2,  tt/2)  and  is  calculated  as 

A  =  sin_1(pz/|p|)  +  ve , 

where  pz  is  the  z-component  in  p  and  ve  is  zero-mean,  white,  measurement  noise.  Care  should 
be  taken  to  avoid  singularities  as  the  target  is  approached  when  \p\  =0.  If  so,  then  ze  =  0.  The 
azimuth  angle  is  measured  on  the  interval  [0,  2i r]  and  is  found  via 

=  atan2 (py,px)  +  va  , 
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where  atan2  is  the  two-argument  arc-tangent  function,  px  and  py  are  the  x  and  y  components  of 
the  relative  position,  p,  and  va  is  zero-mean,  white,  measurement  noise. 

Reference  measurements  can  be  determined  by  evaluating  the  elevation  and  azimuth  angle 
models  at  a  reference  state.  These  reference  measurements  are  then  used  to  define  measurement 
deviations  as  the  difference  between  the  measurement  and  its  reference.  A  first-order  Taylor  series 
expansion  for  the  deviations  in  the  elevation  and  azimuth  angles  yields 

e  <9e  ,  .  da  . 

dze  =  —0 p  +  ve  and  dza  =  - —dp  +  va  , 

op  op 

where  (de/dp)  and  (da /dp)  are  simple  derivatives  of  the  trigonometric  functions  and  5p  is  the  error 
in  the  relative  position  vector  with  respect  to  the  reference.  Depending  on  the  IMU  mechanization, 
this  error  changes  form. 

Strapdown  For  the  non-gimbaled  mechanization,  the  relative  position  vector,  in  inertial  coordi¬ 
nates,  to  the  destination  target  from  the  sensor  is  given  by 

pit/s  =  T}rft-(riMV  +  Tirbs/lMV), 

which  can  be  used  to  determine  the  error  in  the  relative  position  as 

dPt/s  =  _^riMu  +  n[rbs/lMVx]6G  , 

where  56  is  the  attitude  error,  which  is  the  error  of  the  body-fixed  frame. 

Space  Stabilized  With  a  gimbaled,  space  stabilized,  platform,  the  relative  position  vector  of 
the  destination  target  with  respect  to  the  sensor  requires  the  added  transformation  between  the 
platform  and  the  body  frame.  The  result  is  given  by 

Pl,,  =  T}r{  -  (rjBU  +  T;Tyt/lm) . 

Evaluating  a  reference  relative  position  vector,  determining  the  error,  and  performing  a  first-order 
Taylor  series  expansion  yields 

sy,  = 

where  56  is  the  attitude  error,  which  is  the  error  of  the  platform-fixed  frame. 

Data  Processing  In  order  to  process  the  bearing,  it  is  required  to  determine  the  Jacobian  of  the 
measurement,  which  takes  on  the  same  form,  independent  of  the  IMU  mechanization,  such  that  it 
can  be  written  as 


IT  _  [  He,r 

0lx3 

0lx3 

H6)g  0iX24 

Hafi  0iX24j 

where 

Her  =  -(de/dp) 

and 

Har  =  -(da/ dp) 

for  both  of  the  IMU  mechanizations, 

=  (de/dp)fl[rbs/lMVx] 

and 

Ha,r  =  ( da/dp)f/[rb/lMUx ] 

for  the  strapdown  IMU  mechanization,  and 

He0  =  (de/dp)f;[Tbprbs/lMVx] 

and 

H^r  =  (da/dp)f;[Tbprbs/lMV 

for  the  space  stabilized  IMU  mechanization. 

23 

Approved  for  public  release;  distribution  is  unlimited. 


3.7  Linear  Covariance  Analysis 

Recall  from  Section  3.5  that  either  of  the  mechanizations  of  the  IMU  can  be  represented  by  a 
nonlinear  dynamical  system;  provided  that  the  position,  velocity,  and  attitude  are  collected  into  a 
state  vector,  x,  the  general  nonlinear  dynamical  system  representing  the  forward  evolution  of  the 
state  for  either  of  the  IMU  systems  is  given  by 


f  (®fc— 1)  ^ki  t k )  j  (^6) 

where  a and  u>k  represent  the  true  vehicle  non-gravitational  acceleration  and  angular  velocity, 
respectively.  For  the  mechanizations  presented  in  Section  3.4,  /(•)  represents  one  of  Eqs.  (10)  or 
Eqs.  (11).  In  addition,  the  true  vehicle  non-gravitational  acceleration  and  angular  velocity  are 
each  corrupted  by  a  set  of  errors  (4>k  for  the  non-gravitational  acceleration  and  xj)]~  for  the  angular 
velocity)  to  produce  the  measurements  of  the  non-gravitational  acceleration  and  angular  velocity, 
which  is  functionally  represented  by  Eq.  (13).  These  parameters  can  also  evolve  in  time,  as  dictated 
by  Eqs.  (14). 

Given  the  dynamical  systems  models  along  with  a  representation  of  the  statistics  of  the  quanti¬ 
ties  involved,  it  is  possible  to  investigate  the  evolution  of  the  uncertainty  for  the  position,  velocity, 
and  attitude  of  the  vehicle.  As  discussed  in  Section  3.5,  one  method  for  investigating  this  uncer¬ 
tainty  evolution  is  through  the  use  of  Monte  Carlo  simulation.  The  Monte  Carlo  method,  while 
highly  flexible,  requires  significant  computational  resources.  An  alternative  method,  which  has  the 
advantage  of  being  less  computationally  demanding,  is  to  use  a  Kalman  filtering-based  approach. 

The  particular  Kalman  filtering-based  approach  used  here  is  known  as  linear  covariance  analysis. 
In  linear  covariance  analysis,  the  nonlinear  systems  describing  the  evolution  of  the  IMU-driven 
systems  are  linearized  about  a  reference  truth  trajectory.  This  has  the  immense  advantage  of 
decoupling  the  uncertainty  propagation  from  any  state  propagation  method,  which  means  that 
only  the  covariance  needs  to  be  propagated  in  order  to  determine  the  uncertainty  of  the  position, 
velocity,  and  attitude  of  the  vehicle.  In  order  to  incorporate  the  uncertainties  present  in  the 
IMU  error  sources,  an  augmented  state  vector,  ajaug,  is  defined  to  be  the  standard  state  vector, 
x,  concatenated  with  all  of  the  accelerometer  and  gyroscope  error  parameters  except  those  that 
are  characterized  as  white-noise  sequences.  All  of  the  white-noise  error  parameters  are  collected 
together  as  the  vector  wavLg.  With  these  definitions,  a  new  nonlinear  dynamical  system  can  be 
constructed  for  the  augmented  state  vector,  which  is  given  by 

®aug,fc  ,/aug(®aug,fc— 1)  ^aug,fc— l)  •  (17) 

where  a?aug,fc  and  a;aug,fc-i  are  the  (augmented)  state  at  times  tk  and  tk-i  respectively. 

This  dynamical  system  describes  the  evolution  of  the  position,  velocity,  and  attitude  of  the 
vehicle,  as  well  as  the  non-white-noise  error  parameters  associated  with  the  IMU.  The  nonlinear 
equation  in  Eq.  (17)  is  then  linearized  about  some  reference  truth,  xt,  similar  to  how  the  extended 
Kalman  filter  is  linearized  about  the  current  state  estimate.  The  resulting  linearized  dynamical 
system  is  used  to  construct  the  covariance  propagation  equation 

Pk  =  Fk-iPk-iFk-i  +  ,  (18) 

where  Pk  represents  the  covariance  at  time  tk-  In  Eq.  (18),  F^- 1  is  the  tangent  linear  dynamics 
evaluated  along  the  reference  state  and  nominal  IMU  model  parameters,  Mk-i  is  the  Jacobian  of  the 
state  with  respect  to  the  white-noise  sequence  parameters,  and  Qk-i  is  the  process  noise  covariance 
matrix  that  represents  the  statistics  of  the  white-noise  sequence.  This  covariance,  represented  by 
Eq.  (18),  represents  the  covariance  of  the  full  augmented  state  vector;  to  investigate  and  analyze 
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the  uncertainties  characterizing  the  position,  velocity,  and  attitude  of  the  vehicle,  the  corresponding 
elements  can  be  extracted. 

The  linear  covariance  approach  can  also  be  used  to  investigate  the  effects  of  acquiring  external 
sensor  data,  such  as  that  from  a  GPS  sensor,  an  altimeter,  a  range  sensor,  or  a  magnetometer, 
among  others.  External  data  is  generally  employed  to  reduce  the  uncertainty  that  is  present  in  the 
prediction  of  the  translational  and  rotational  states  of  a  vehicle.  By  combining  external  data  with 
IMU  data,  lower  uncertainties  can  be  obtained;  however,  this  means  relying  on  data  regarding  the 
external  environment,  which  can  lead  to  sensitivities  in  data  quality  and  reliability.  The  external 
data  are,  in  general,  taken  to  be  of  the  form 

%k  —  ^-(®aug,fc)  T  Vk  , 

where  vk  represents  the  measurement  noise  of  the  sensor  at  time  tk,  which  is  taken  to  be  zero  mean 
with  covariance  Rk.  When  presented  with  this  external  measurement  data,  the  linear  covariance 
method  updates  the  covariance  through  a  blending  of  the  a  priori  state  information  and  the  newly 
acquired  external  measurement  data  via 

P+  =  (I-  KkHk)  Ph7  ( I  -  KkHk)T  +  KRkKT  , 

where  Kk  is  the  Kalman  gain,  which  is  given  by 

Kk  =  Pk  Hj.  ( HkP~Hl  +  Rkyl  . 

Whereas  the  nonlinear  dynamics  function,  /aUg(-)>  and  the  tangent  linear  dynamics  matrix,  Fk~i, 
are  key  in  the  propagation  stage  of  the  linear  covariance  method,  the  nonlinear  measurement 
function,  h(-)  and  its  associated  tangent  linear  matrix,  Hk ,  are  critical  for  the  update  stage  of  the 
linear  covariance  method. 

4.0  RESULTS  AND  DISCUSSION 

4.1  Trajectory  Simulation 

In  order  to  provide  the  trajectory  input  to  SAIMUN  that  is  outlined  in  Section  3.1,  a  ground 
vehicle  was  used  to  log  acceleration  and  angular  velocity  data  using  a  VectorNav  VN-1001  IMU. 
The  VN-100  logged  data  at  40  Hz  during  a  trajectory  that  took  approximately  41  minutes  to 
complete.  A  smoother  was  implemented  to  decrease  the  effect  of  noise  on  the  recorded  acceleration 
and  angular  velocity  data.  The  smoothed  acceleration  and  angular  velocity  data  were  dead-reckoned 
from  a  known  starting  condition  to  provide  the  true  trajectory  for  SAIMUN  to  use.  The  dead- 
reckoned  trajectory  is  not  perfectly  representative  of  the  trajectory  followed  by  the  vehicle;  however, 
for  the  purposes  of  applying  SAIMUN  and  analyzing  IMU-based  navigation  performance,  the  dead- 
reckoned  trajectory  is  taken  to  be  the  true  trajectory. 

4.2  IMU  Selection 

To  demonstrate  the  utility  of  SAIMUN,  two  inertial  measurement  units  are  compared:  the  Lord 
Microstrain  3DM-GX3-152  and  the  Analog  Devices  ADIS164883.  Both  IMUs  utilize  MEMS  sensors 

1  http://www.vectornav.com/Downloads/Support/PB-12-0002.pdf,  accessed  04/21/2014 
2 http: //files. microstrain.com/3DM-GX3- 15- Inertial- Measurement-Unit. pdf,  accessed  04/21/2014 
3http://www. analog.com/static/imported-files/data_sheets/ ADIS  16488. pdf,  accessed  04/21/2014 
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to  provide  triaxial  accelerometers,  gyroscopes,  magnetometers,  and  a  pressure  sensor.  The  IMUs 
are  implemented  to  provide  dead-reckoning  navigation  of  a  ground  vehicle  trajectory. 

For  this  analysis,  only  the  accelerometer  and  gyroscopic  information  is  of  interest.  Table  1  shows 
the  pertinent  specifications  for  both  of  the  chosen  IMUs.  These  IMUs  are  interesting  to  compare 
due  to  the  differences  between  the  random  walk  and  bias  instability  specifications.  Table  1  shows 
that  the  ADIS  16488  has  a  higher  random  walk  specification  and  thus  more  white  noise  present 
in  the  signal.  The  3DM-GX3-15,  on  the  other  hand,  has  a  higher  bias  instability  specification. 
The  higher  bias  instability  and  lower  random  walk  specifications  of  the  3DM-GX3-15  compared  to 
the  ADIS16488  imply  that  the  position  and  attitude  uncertainty  associated  with  the  3DM-GX3-15 
will  grow  more  slowly  than  the  ADIS  16488  initially.  Eventually,  the  position,  velocity,  and  atti¬ 
tude  uncertainty  associated  with  the  3DM-GX3-15  will  grow  larger  and  surpass  their  counterparts 
associated  with  the  ADIS  16488  as  time  increases. 

Table  1.  IMU  Specifications  for  the  ADIS16488  and  3DM-GX3-15 


Specification 

3DM-GX3-15 

ADIS16488 

Analog-to-Digital  Converter  Bitrate 

16 

32 

Range 

±156.96  f 

±176.58  f 

c n 

pH 

CD 

Scale  Factor  Error 

Z4(— 500,500)  ppm 

TV ( — 5000,  5000)  ppm 

CD 

a 

Axes  Nonorthogonality  Error 

Z4(— 103, 103)  arcsec 

TV ( — 126, 126)  arcsec 

o 

Ph 

JD 

Axes  Misalignment  Error 

Z4(— 103, 103)  arcsec 

U{— 3600,  3600)  arcsec 

"3 

CD 

CD 

Velocity  Random  Walk 

0.0007848^2 

Vhr 

0.029 

Vhr 

<1 

Bias  Instability 

0.0003924  f 

0.000981  f 

Bias  Instability  Time 

100  s 

100  s 

Startup  Bias 

«(-0.0196, 0.0196)  ™ 

TV ( — 0.016,  0.016)  f 

Analog-to-Digital  Converter  Bitrate 

16 

32 

Range 

±600 

±450 

CO 

Scale  Factor  Error 

Z4(— 500,  500)  ppm 

TV( — 100,  100)  ppm 

CD 

a 

o 

Axes  Nonorthogonality  Error 

U{— 103, 103)  arcsec 

TV ( — 180, 180)  arcsec 

CD 

CO 

O 

Axes  Misalignment  Error 

TV ( — 103, 103)  arcsec 

U[— 3600,  3600)  arcsec 

o 

Angular  Random  Walk 

0.03  % 

Vhr 

0.3 

Vhr 

Bias  Instability 

18  If 

6.25^f 

Bias  Instability  Time 

100  s 

100  s 

Startup  Bias 

Z4(— 0.25,0.25)  ^ 

TV ( — 0.2,  0.2) 

4.3  Performance  Metrics 

In  order  to  assess  the  performance  of  an  IMU  for  a  given  trajectory  the  standard  deviation  in  the 
position,  velocity,  and  attitude  is  utilized  as  the  performance  metric.  Because  the  ADIS16488  and 
the  3DM-GX3-15  are  both  relatively  low-end  IMUs,  two  cases  are  considered:  Case  1  analyzes  only 
the  first  500  seconds  of  the  entire  41  minute  trajectory  while  Case  2  analyzes  the  entire  41  minute 
trajectory. 
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Figure  12.  Position  Standard  Deviation  as  a  Function  of  Time  for  the  ADIS16488  and 
the  3DM-GX3-15  IMUs:  Case  1 

4.3.1  Standard  Deviation  in  Position  and  Attitude.  In  order  to  quantify  uncertainty  in 
the  position,  velocity,  and  attitude  states,  the  standard  deviation,  or  square  root  of  the  second- 
central  moment,  is  calculated  and  plotted  as  a  function  time  based  on  the  Monte  Carlo  samples 
or  the  linear  covariance  estimate.  The  position  and  attitude  standard  deviations  computed  for  the 
ADIS16488  and  3DM-GX3-15  IMUs  can  be  found  in  Figs.  (12)  and  (13)  respectively  for  the  Case 
1  trajectory  with  velocity  omitted  for  brevity.  Due  to  the  similarity  between  the  strapdown  and 
space  stabilized  metrics,  only  the  strapdown  results  are  presented.  As  was  stated  in  Section  4.2, 
the  ADIS  under-performs  against  the  3DM  at  the  beginning  of  the  run,  but  this  switches  within 
the  Case  1  bounds.  Through  the  rest  of  the  run,  the  ADIS  continues  to  outperform  the  3DM.  The 
Case  2  plots  are  omitted  for  the  IMU  comparison  for  conciseness;  after  the  500  second  mark,  the 
two  IMUs  continue  to  diverge. 

To  test  the  validity  of  linear  covariance  analysis,  the  standard  deviation  of  the  position,  velocity, 
and  attitude  for  the  strapdown  mechanization  are  plotted  against  the  same  metric  computed  via 
Monte  Carlo  simulation.  The  results  for  Case  2  are  shown  in  Figs.  (14)— (16).  With  the  exception 
of  the  z  component  of  the  velocity  in  Fig.  (15),  the  linear  covariance  and  Monte  Carlo  analyses 
match  very  closely.  The  standard  deviations  observed  from  Case  1  are  so  nearly  identical,  the  plots 
were  not  included. 

Similarly  to  the  strapdown  mechanization,  the  space  stabilized  results  from  linear  covariance 
and  Monte  Carlo  match  quite  closely  with  the  exception  of  the  velocity  z  component.  When 
scrutinized,  the  deviation  of  linear  covariance  from  Monte  Carlo  is  smaller  when  the  space  stabilized 
mechanization  is  implemented.  This  is  to  be  expected  due  to  the  nature  of  the  dynamics.  Newtonian 
mechanics  produce  the  simplest  result  when  the  acceleration  components  occur  naturally  in  the 
inertial  frame  which  is  the  intent  of  a  space  stabilized  mechanization.  Conversely,  the  accelerations 
obtained  via  a  strapdown  mechanization  must  be  converted  from  the  body  to  the  estimated  inertial 
frame. 
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Figure  13.  Attitude  Standard  Deviation  as  a  Function  of  Time  for  the  ADIS16488  and 
the  3DM-GX3-15  IMUs:  Case  1 
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Figure  14.  Position  Standard  Deviation  for  the  ADIS  16488  in  a  Strapdown  Mecha¬ 
nization:  Case  2 
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Figure  15.  Velocity  Standard  Deviation  for  the  ADIS16488  in  a  Strapdown  Mecha¬ 
nization:  Case  2 
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Figure  16.  Attitude  Standard  Deviation  for  the  ADIS16488  in  a  Strapdown  Mecha¬ 
nization:  Case  2 
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Figure  17.  Position  Standard  Deviation  for  the  ADIS16488  in  a  Space  Stabilized 
Mechanization:  Case  2 
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Figure  18.  Velocity  Standard  Deviation  for  the  ADIS16488  in  a  Space  Stabilized 
Mechanization:  Case  2 
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Figure  19.  Attitude  Standard  Deviation  for  the  ADIS16488  in  a  Space  Stabilized 
Mechanization:  Case  2 


Table  2.  External  Aiding  Configuration 


Sensor 

Start  Time  (s) 

Finish  Time  (s) 

Meas.  Spacing  (s) 

Sensor  Noise 

Range/Range- Rate 

0 

150 

10 

0.1  m-m/s 

Stellar  Line-of-Sight 

400 

600 

10 

20  arcsec 

Bearing  Angles 

800 

1000 

0.1 

20  arcsec 

4.3.2  Standard  Deviation  with  External  Aiding.  The  effects  of  the  three  external  sensors 
(range/range-rate,  stellar  line-of-sight,  and  bearing  angles)  on  the  standard  deviation  growth  for 
the  space  stabilized  mechanization  are  presented  in  Figs.  (20)-(22).  The  sensors  are  run  according 
to  the  values  given  in  Table  2.  The  figures  show  the  ADIS  with  varying  levels  of  external  aiding: 
unaided,  range/range-rate  and  stellar  line-of-sight,  and  stellar-line-of-sight  and  bearing  angles.  It 
should  be  noted  that  the  noise  associated  with  each  of  the  sensors  is  assumed  to  be  zero-mean 
white  noise. 

Based  on  the  standard  deviation  growth  pattern,  although  the  late  stage  bearing  angles  sig¬ 
nificantly  reduce  translational  uncertainty,  initial  range/range-rate  measurements  are  much  more 
beneficial. 

As  a  demonstration  of  the  flexibility  of  the  tool,  Figs.  (23)-(25)  depict  a  hypothetical  case  study 
attempting  to  determine  a  combination  of  sensors  that  produces  standard  deviation  growth  similar 
to  a  “tactical”  grade  IMU.  Shown  is  an  unaided  fictional  tactical  IMU  (achieved  by  reducing  error 
sources  in  the  ADIS  IMU),  an  unaided  ADIS,  and  an  aided  ADIS.  Based  on  this  analysis,  a  user 
could  determine  that  the  use  of  the  three  presented  external  sensors  could  ensure  the  uncertainty 
growth  associated  with  an  ADIS  IMU  would  be  comparable  to  the  tactical  IMU.  A  user  could  then 
confidently  assume  that  the  purchase  of  a  costly  high-grade  IMU  would  be  unnecessary  provided 
that  the  external  sensor  data  is  available  via  the  additional  sensors.  These  results  are  easily  obtain¬ 
able  with  minimal  computational  effort.  Had  a  case  study  of  this  sort  been  attempted  in  a  Monte 
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Figure  20.  Position  Standard  Deviation  from  the  ADIS16488  in  a  space  stabilized 
mechanization  with  and  without  external  aiding.  EA1  utilizes  range/range-rate  and 
stellar  line-of-sight  measurements,  while  EA2  utilizes  stellar  line-of-sight  and  bearing 
angles. 
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Figure  21.  Velocity  Standard  Deviation  from  the  ADIS16488  in  a  space  stabilized 
mechanization  with  and  without  external  aiding.  EA1  utilizes  range/range-rate  and 
stellar  line-of-sight  measurements,  while  EA2  utilizes  stellar  line-of-sight  and  bearing 
angles. 
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Figure  22.  Attitude  Standard  Deviation  from  the  ADIS16488  in  a  space  stabilized 
mechanization  with  and  without  external  aiding.  EA1  utilizes  range/range-rate  and 
stellar  line-of-sight  measurements,  while  EA2  utilizes  stellar  line-of-sight  and  bearing 
angles. 


Carlo  simulation  (with  external  aiding  utilities),  this  process  would  have  required  a  much  larger 
computational  burden. 
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Figure  23.  Position  standard  deviations  resulting  from  an  unaided  tactical  IMU  and 
the  ADIS16488  IMU  with  and  without  external  aiding  in  a  space  stabilized  mecha¬ 
nization.  All  three  sensors  listed  in  Table  2  are  utilized. 
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Figure  24.  Velocity  standard  deviations  resulting  from  an  unaided  tactical  IMU  and 
the  ADIS16488  IMU  with  and  without  external  aiding  in  a  space  stabilized  mecha¬ 
nization.  All  three  sensors  listed  in  Table  2  are  utilized. 


34 

Approved  for  public  release;  distribution  is  unlimited. 


I  I  I  r  ...............  I  I  -  I  I 

0  100  200  300  400  500  600  700  800  900  1,000 


Time  [s] 


Figure  25.  Attitude  standard  deviations  resulting  from  an  unaided  tactical  IMU  and 
the  ADIS16488  IMU  with  and  without  external  aiding  in  a  space  stabilized  mecha¬ 
nization.  All  three  sensors  listed  in  Table  2  are  utilized. 
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5.0  CONCLUSIONS 


Before  choosing  an  IMU  and  external  sensors  to  be  implemented  in  a  navigation  system,  SAIMUN 
can  be  used  to  analyze  the  performance  of  the  full  sensor  suite  along  a  given  trajectory  or  trajec¬ 
tories.  When  comparing  the  Lord  MicroStrain  3DM-GX3-15  and  the  Analog  Devices  ADIS16488 
IMUs  in  the  navigation  of  an  example  vehicle  trajectory,  the  lower  random  walk  specification  of 
the  3DM-GX3-15  causes  its  associated  uncertainty  to  grow  at  a  slower  rate  during  the  beginning 
of  the  trajectory;  however,  the  higher  bias  instability  of  the  3DM-GX3-15  causes  its  associated 
uncertainties  to  grow  faster  and  overtake  the  uncertainties  associated  with  the  ADIS16488.  For 
these  IMUs,  425-500  seconds  is  observed  to  be  the  approximate  window  in  which  the  performance 
of  the  IMUs  is  most  similar,  as  quantified  by  the  standard  deviations.  Before  this  window,  the 
3DM-GX3-15  provides  a  more  accurate  navigation  solution.  After  this  window,  the  ADIS16488 
provides  a  more  accurate  navigation  solution.  SAIMUN  allows  this  trade  study  to  be  performed 
before  purchasing  or  testing  any  hardware  and  will  allow  for  larger  trade  spaces  to  be  investigated 
virtually,  obtaining,  building,  or  testing  any  physical  IMUs. 

SAIMUN  enables  both  gimbal-  and  strap  down-based  navigation  systems  to  be  analyzed  virtually 
in  the  developed  simulation  environment.  An  IMU  error  parameterization  and  a  trajectory  are 
defined  in  the  simulation  environment  along  with  a  variety  of  external  sensors,  including  their  own 
error  parameterizations.  To  provide  validity  of  the  unaided  linear  covariance  analysis,  both  linear 
covariance  and  Monte  Carlo  analyses  can  be  run  in  the  same  simulation.  With  linear  covariance 
analysis  verified  via  Monte  Carlo  analysis,  the  external  aiding  sensors  can  be  added  to  the  simulation 
and  Monte  Carlo  can  be  removed. 

Currently,  SAIMUN  is  capable  of  rendering  navigation  solutions  for  gimbal-  and  strapdown- 
based  navigation  systems,  while  considering  IMU  error  sources  due  to  startup  bias,  bias  instabil¬ 
ity,  thermo-mechanical  zero-mean  white  noise,  scale  factor  errors,  axes  misalignment  errors,  axes 
nonorthogonality  errors,  and  quantization  effects  caused  by  analog-to-digital  conversion.  Where 
appropriate,  SAIMUN  allows  for  a  selection  from  multiple  probabilistic  representations  of  the  IMU 
error  sources,  such  as  Dirac,  Gaussian,  and  uniform  distributions.  The  structure  of  SAIMUN  sup¬ 
ports  future  extensions,  such  as  additional  IMU  mechanizations,  IMU  parameter  errors,  error  source 
distributions,  and  sensors.  On  the  analysis  side  of  SAIMUN,  this  is  also  true;  that  is,  analyses  ex¬ 
tending  beyond  those  presented  in  this  report  can  be  incorporated.  The  simulation  and  analysis 
environment  is  designed  to  be  modular  such  that  these  modifications  are  straightforward  to  imple¬ 
ment  within  SAIMUN. 

Multiple  possibilities  for  external  aiding  are  incorporated  into  SAIMUN,  such  as  range/range- 
rate  measurements,  stellar  line-of-sight  measurements,  and  bearing  angles  measurements.  SAIMUN  fa¬ 
cilitates  analysis  that  combines  both  IMU  and  external  aiding  data  to  investigate  the  possible  nav¬ 
igation  accuracy  of  aided  or  unaided  navigation  solutions.  In  this  way,  SAIMUN  can  be  used  to 
determine  if  a  lower  cost  IMU  can  be  paired  with  a  set  of  external  sensors  to  achieve  similar  or 
better  performance  as  a  higher  cost  IMU.  It  is  shown  that  combining  all  three  external  sensors  with 
a  MEMS  IMU  can,  for  the  trajectory  investigated,  lead  to  similar  navigation  accuracies  as  would 
be  achieved  with  a  fictitious  tactical  grade  IMU. 
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LIST  OF  SYMBOLS,  ABBREVIATIONS,  AND  ACRONYMS 


AIRS  Advanced  Inertial  Reference  Sphere 
CM  Center  of  Mass 
EKF  Extended  Kalman  Filter 
FPA  Focal  Plane  Angles 
IMU  Inertial  Measurement  Unit 


MEMS  MicroElectricalMechanical  Systems 
ng  non-gravitational 
pdf  Probability  Density  Function 
QUEST  Quaternion  Estimation 
SAIMUN  Simulation  and  Analysis  of  IMU-based  Navigation 
VRW  Velocity  Random  Walk 
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